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SUMMARY 


Because of the complexity of the helicopter flow field, a zonal method 
of analysis of computational aerodynamics is required. The objective of 
this research is to demonstrate the feasibility of a new potential/viscous 
coupling procedure for reducing computational effort while maintaining solu- 
tion accuracy. 

A closed- loop overlapping velocity- coupling procedure has been devel- 
oped with the unique feature that the potential flow singularity strengths 
are obtained from the Navier- Stokes solution at an inner fluid boundary. 
This procedure has been utilized to combine a two-dimensional potential flow 
code and a two-dimensional Navier-Stokes code, ARC2D, developed by Pulliam 
and Steger at NASA Ames. 

The resultant fully coupled code, ZAP2D (Zonal Aerodynamics Program- 
2D) , has been used to compute the flow past a NACA 0012 airfoil. For this 
case the zonal method has shown that the grid domain size can be reduced to 
0.25 chord lengths with no loss in accuracy. The numerical behavior of the 
ARC2D solution is also greatly improved due to the reduced domain size. 
Finally, the computation time is reduced by approximately a factor of 10 for 
given accuracy. 

Additional reductions in domain size and accelerated convergence 
methods should decrease the required ZAP2D computation time to about a fac- 
tor of 20 compared with the ARC2D simulation. This level of computational 
speed achievable with ZAP2D would allow for the development of a fully in- 
teractive program that could be used for preliminary airfoil design. More 
importantly, perhaps, similar grid density reductions in three dimensions 
would result in significant decreases in memory requirements and potential 
computation time reductions for such a zonal approach of at least one order 
of magnitude. 

This improvement in computational effort should then provide the means 
for more accurate viscous and compressible simulations of the complex flow 
fields associated with rotary wing aircraft and should translate into im- 
proved analysis of such important performance parameters as helicopter fuse- 
lage drag and rotor airloads. 
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1.0 INTRODUCTION 


The flow fields associated with rotary wing aircraft, including the 
tilt-rotor helicopter, are extremely complicated and remain a challenge for 
applied computational aerodynamics. While linear methods or pure potential 
flow models of nonlinear wake shape effects have become commonly used tech- 
niques for helicopter aerodynamics predictions and analyses, the nonlinear 
methods that incorporate more correct mathematical models of the actual flow 
physics- -such as compressibility and viscosity- -are being developed and ap- 
plied by research scientists. Much progress has been made in solving the 
full potential, Euler, and various forms of the approximate Navier- Stokes 
equations; however, because of the wide variation of physical scales of 
characteristic fluid phenomena from boundary layer thickness and vortex core 
radius to rotor diameter, a zonal method of analysis is generally projected 
for the entire helicopter flowfield calculation. This is necessary not only 
to achieve a closed calculation on the next generation of computers, but 
also to obtain the required numerical accuracy over such large regions that 
must include fine meshes to resolve the viscous layers as well as large flow 
gradients. Consequently, full closed-loop coupling methods need to be de- 
veloped for these zonal calculations. 

In the Phase I work described in this report, a new general procedure 
for coupling potential and viscous flow calculation schemes has been devel- 
oped and demonstrated for a two-dimensional airfoil simulation. The coup- 
ling method is based on the premise that any computational method can only 
produce valid results within the approximations of the physical model em- 
ployed in its construction. Therefore, the interfacing boundary surface be- 
tween the potential and viscous flow code domains must be very nearly poten- 
tial in order to achieve a successful (accurate) result. Specifically, dur- 
ing Phase I, an AMI (Analytical Methods, Inc.) two-dimensional, potential 
flow panel method, program P0T2D , has been coupled with a NASA Ames two-di- 
mensional, thin-layer, Navier-Stokes code, ARC2D (1). An overlapping veloc- 
ity coupling procedure is utilized with the unique feature that (beyond the 
first iteration) the potential flow surface panel strengths (^ and d<f>/d n) 
are obtained directly from the Navier-Stokes solution at a smooth inner 
fluid boundary. These fluid surface panel values are then used to compute 
the outer velocity boundary conditions for the Navier-Stokes calculation. 
The iteration between the potential flow solution and the Navier-Stokes so- 
lution continues in a closed loop until flow convergence or allowable itera- 
tion limits are obtained. The resultant fully coupled code is called ZAP2D 
for £onal Aerodynamics program (2D). 

The primary objective of the current work is then to demonstrate the 
feasibility of the coupling procedure installed in ZAP2D for reducing the 
computational effort, i.e., grid, required by ARC2D for a given accuracy. 
In the following sections the methodology required for the coupled program 
is first described briefly and then the detailed analysis of the zonal simu- 
lation of the flow past an airfoil is compared with the Navier-Stokes alone 
simulation as well as with experimental data. 
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2.0 


METHODOLOGY 


2 . 1 Zonal Method- -ZAP2D 

A new two-dimensional CFD method, ZAP2D (Zonal Aerodynamics Program), 
has been developed that combines potential and viscous flow calculation 
schemes. The actual coupling between the methods is based on the premise 
that the coupling boundary should be located in the computed flow field 
where the approximations inherent in both methods are valid; consequently, 
the interfacing boundary surface between the potential and viscous flow re- 
gions must be very nearly potential in order to achieve an accurate result. 

The general coupling concept for the zonal method is illustrated in 
Figures 1 and 2. Figure 1 shows the zonal representation of the physical 
flow field while the flow chart of the iteration scheme is shown in Figure 

2. An inner boundary, which is measured by a radius, say, R. , is con- 

structed to enclose a generally shaped body or airfoil and the iluid region 
near the body where viscous effects dominate. Although the body might be a 
complex shape, the inner boundary surface, S., would be a simple, smooth 
geometry. An additional outer surface, S , whicft is measured by a radius, 
say, R q , also of simple shape, forms the computational domain of the viscous 
calculation. Of course, the potential flow is solved by the panel method 
beyond S q to an infinite distance from the body surface. 

The iterative coupling of the proposed method proceeds as follows (see 
Figure 2) . 

1. The calculation procedure is initialized by computing the potential 
flow associated with the body motion. An integral panel method is used 
for this computation. In this first step only the actual body surface 
is represented by discrete panels and the surface potential equation is 
solved directly by enforcing the Neumann boundary condition. 

2. The velocity field, V , on the outer boundary, S , is computed from the 
known potential flow surface singularity solution. 

3. This velocity field, V , serves as the outer boundary values for the 
Navier- Stokes calculation. Once the Navier- Stokes solution for the 
flow inside S q is obtained, the velocity field due to the viscous solu- 
tion, say V. on S^, which is just outside the region of viscous ef- 
fects, is also known. 

4. This inner boundary velocity, V. , can therefore be used to generate the 

known corresponding values of the velocity potential, <f>, and its normal 
derivative to S., (d<f>/d n) . These represent the panel doublet and 

source strengths on required for a new potential flow calculation; 
consequently, all further potential flow calculations require only the 
paneling of the smooth inner boundary, S., and not the actual body sur- 
face. Therefore, further calculations by the panel method of the outer 
boundary velocity field, V , should be very accurate and any possible 
contamination of the Navier-stokes simulation by potential flow numeri- 
cal errors is minimized. 
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Figure 1. Zonal Representation of the Physical Flow Field 
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Figure 2. Flow Chart of the Zonal Iteration Scheme. 
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This iteration loop, Steps 2-4, is repeated until the inner boundary 
velocity, V^, is converged as measured by converged integrated forces 
and moments. 

In order to validate such a coupling procedure it is important to se- 
lect robust and computationally accurate potential and viscous simulation 
methods. P0T2D was selected as the potential flow module since it was de- 
veloped by AMI and formed the basis of the three-dimensional code, VSAERO 
(2), which has been widely used and validated. ARC2D was selected for the 
viscous flow module for much the same reasons- -it is very widely used in the 
computation of airfoil aerodynamics and has also been extended to the three- 
dimensional case, ARC3D. 

These two codes, P0T2D and ARC2D, are briefly described below with em- 
phasis given to the code modifications and accuracy testing required to ac- 
complish the coupling procedure. In particular, the validity of the poten- 
tial flow simulation for the open inner boundary, S., is of importance and 
is described in some detail. 1 


2 . 2 Potential Flow Method- -P0T2D 

As is well known, the potential describing an inviscid, irrotational , 
incompressible flow can be expressed by its boundary values (3) . Specifi- 
cally, an integral solution of Laplace's equation for the velocity potential 
can be derived by application of Green's theorem. The resultant boundary 
integral equation for the surface velocity potential is composed of two in- 
tegrals: 

1. a body surface integral of a doublet kernel function of strength, <f> t 
and a source kernel function of strength, d<f>/d n, and 

2. a thin wake surface integral of a doublet kernel function of 
strength, A<^. 

Of course, the wake strengths, A <f> , are related to the body surface poten- 
tial, through the wake auxiliary Kutta condition and the body source 
strengths, d<f>/d n, are known by the tangent flow boundary condition on the 
body. Hence the integral equation for <f> can be solved directly as suggested 
by Morino (4). In P0T2D, the integral equation is discretized and trans- 
formed to an algebraic equation by representing the boundaries with flat 
panels and approximating the local panel source and doublet strengths with 
constant values just as in VSAERO. 

Once the surface potential is computed, surface velocities are obtained 
by numerical differentiation and forces and moments are found by pressure 
integration. Further, off -body flow velocities are computed by integration 
of surface velocity source and doublet influence functions. Finally, in the 
current code, the Karman-Tsien compressibility correction (5) based upon 
free stream Mach number and local velocity magnitude can also be utilized. 
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In the continued development of P0T2D for the zonal method, surface 
solution accuracy, off -body velocity scan accuracy, and the accurate simula- 
tion of the open surface required for the inner boundary coupling were care- 
fully analyzed. A summary of the specific test cases used and their impli- 
cations are presented in Table 1. The results of each of these cases are 
discussed briefly below. 

In case 1, the accuracy of the P0T2D algorithm was tested by comparison 
with the exact solution for a circular cylinder with a dimensionless circu- 
lation equal to approximately 2.2. As shown in Figure 3, the comparison of 
the computed and exact surface pressure coefficient is remarkably accurate 
for a panel model consisting of only 80 upper surface elements. The com- 
puted off-body velocity magnitude is compared to the exact solution in 
Figure 4. The closest field velocity calculation is located approximately 
10% of a diameter (or chord length) away from the cylinder surface. This 
distance would represent a typical minimum location for the inner boundary 
coupling surface. In this calculation, the panel number was also increased 
to 240 equally spaced panels in order to simulate the relative ratio of R. 
to the maximum local panel size that would be expected for the airfoil 
calculations. As indicated in Figure 4, the maximum deviation in velocity 
magnitude from exact value at this location is on the order of one percent, 
which should be suitable for the zonal method. 

A final example of the accuracy of P0T2D for simulating the flow past a 
closed surface is shown in Figure 5. Here, the computed and experimental 
(6) surface pressures are compared for the NACA 0012 airfoil at 4.966 de- 
grees angle of attack (Case 2). As expected, the potential flow solution 
from P0T2D is very accurate over most of the airfoil, with the exception of 
areas near the leading and trailing edges, where viscous effects are impor- 
tant . 

Since the P0T2D simulation at the fluid inner boundary coupling sur- 
face, S^, is driven by the ARC2D solution which is computed only a finite 
distance downstream, the potential flow simulation on this open surface must 
be extended to infinity. A wake model which implies a branch cut to infini- 
ty - was installed in P0T2D which cancels the vortex left at the upper and 
lower trailing edges due to the truncation of the potential doublet panels. 
A simple Rankine oval was utilized to verify the technique (Case 3). As 
presented in Figure 6, the deviation from the exact surface pressure coeffi- 
cient is negligible over the entire surface and, in fact, is only evident as 
the trailing edge of the body is approached. This small error near the open 
end has recently been found to be due to the truncation of the perturbation 
source strength at the rear of the surface. While the influence of a linear 
source panel has been derived and added to the code to carry the source 
strength off onto the wake in order to recover free -stream conditions, this 
contribution has not been properly verified. Consequently, in all calcula- 
tions presented here, this term is neglected, and the accuracy demonstrated 
here is assumed adequate for the proof-of-concept. Still, future improve- 
ments should include the proper decay of the surface source term to ensure 
an accurate solution which would be independent of the location of the inner 
boundary . 
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DATA CORRELATION 

CASE DESCRIPTION PURPOSE EXACT OR EXP. ALPHA MACH RESULTS IN 

NUMERICAL NO. FIG. NO. 
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Table 1. Summary of Test Cases used in the Potential Flow Validation. 
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Figure 4. Calculated Field Velocity Error (%) from the Exact Potential Flow Solution for Two- 
Dimensional Circular Cylinder (Case 1, Table 1). 
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Figure 5. Comparison of Calculated and Experimental Cp Distribution for 

a NACA 0C12 Airfoil at Mach 0.3, Alpha 4.966° (Case 2, Table 1). 
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Figure 6. Comparison of Calculated and Exact Potential Flow Solution of Cp 
Distribution for Rankine Body (Case 3, Table 1). 
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To verify the mechanics of the ZAP2D coupling procedure, a rather in- 
teresting test case was devised (Case 4) . This case involved a pseudo- 
coupling procedure which was purely potential flow in nature. The initial 
part of this case involved the surface solution for the NACA 0012 airfoil at 
5° angle of attack and a field velocity computation on a surface much like 
the inner boundary (S^) shown in Figure 1. The field velocity components 
were combined to form the normal gradient in potential and used as boundary 
conditions for the second part of the validation case . This "inner bound- 
ary" solution, complete with d$/di\ specification, onset flow, and simple 
wake extension was then compared to the initial field velocity calculation 
in terms of the velocity components. As demonstrated by Figures 7(a) and 
(b) , the part two inner boundary solution compares extremely well with the 
initial field solution for both V and V , respectively. The minor devia- 
tions near the peak velocities and near Z the wake region may be explained by 
the finite representation of the airfoil surface as this relates to the 
field velocity computation and the premature truncation of the trailing-edge 
source term previously mentioned. 

The comparison between the field velocity calculation in the presence 
of the closed surface of part one (ITER-0) and the open surface (ITER=1) of 
part two is presented in Figures 8(a) and (b) for the two components of 
velocity. As indicated, the field velocity calculation generates nearly 
identical results regardless of whether the surface is open or closed. 

As demonstrated by this battery of test cases, the accuracy of the po- 
tential flow solution generated by P0T2D for both closed and open surfaces 
should minimize any induced numerical errors in the ARC2D boundary condi- 
tions . 


2 . 3 Navier- Stokes Method- -ARC2D 

As noted earlier, ARC2D, originally developed by Pulliam and Steger (1) 
at NASA Ames, has been widely used in the computation of airfoil aerodynam- 
ics. The code has been continually improved by Pulliam and colleagues (7), 
(8) and (9), and the details of the theory are well documented in these ref- 
erences. Even though this code is written based on the "thin layer Navier- 
Stokes approximation" it is known to be adequate for separation problems 
with subsonic or transonic free streams when the flow field associated with 
separation is convection-dominated. 

For turbulent flows, the well-known Baldwin-Lomax algebraic model (10) 
is used for turbulence closure. This turbulence model has been used in com- 
puting solutions for a wide variety of flow conditions and has been found to 
be acceptably accurate. 

The numerical algorithm is an implicit approximate factorization finite 
difference scheme which can be either first- or second-order accurate in 
time . 
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(b) Vertical Components 

Figure 7. Comparison of Velocity Components for the Field Calculation 
and the "Equivalent" Surface Solution (Case 4, Table 1). 
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Figure 8. Comparison of Velocity Components for the Field Velocity Calculation 
on the Outer Boundary in the Presence of Both the Airfoil (ITER=0) 
and the "Equivalent" Surface (ITER=1). 
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Over the years of development of ARC2D, unique features have been added 
to the code to improve its capability to accurately predict aerodynamic 
characteristics for airfoil calculations. 

One of these features is a far- field boundary treatment using the 
Riemann invariants concept and the other is an implementation of a simple 
viscous- inviscid coupling. In the far- field boundary condition, the outer 
boundary velocities are determined as a function of both the outer boundary 
velocities and the velocities extrapolated from the inner points. In this 
calculation the normal velocity component at the boundary is determined as 
an average of the Riemann invariants. The tangential velocity component is 
obtained from either the outer boundary velocities at an inflow condition or 
the extrapolated velocities at an outflow condition. Additionally, a simple 
point vortex correction model has been installed to reduce the solution de- 
pendency on the size of the outer boundary. Here, circulation on the far- 
field boundary is set based on the ARC2D computed C.. Using this velocity 
correction together with a compressibility correction, modified outer veloc- 
ities are constructed. It is noted that this correction yields accurate re- 
sults for domains larger than 4 chords in length (11). 

In ZAP2D, a parallel procedure is employed whereby the P0T2D computed 
velocities (currently without compressibility correction) are used as far- 
field boundary conditions. 
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3.0 COMPUTED RESULTS 


The effectiveness of the coupling procedure installed in ZAP2D is meas- 
ured by demonstrating accurate and converged solutions for decreasing domain 
size. In this work, the NACA 0012 airfoil was selected for the validation 
cases because of the experimental data base (6) . 

To complete the schedule of computer simulations required in this 
study, ARC2D and ZAP2D were ported to the Ardent TITAN and also to the IRIS 
4D/20G and IRIS 4D/80GT workstations. In addition, graphics routines were 
also installed on these workstations in order to analyze the results. Fur- 
thermore, a new grid generation routine (GRD2D) was written for these ma- 
chines that can be eventually embedded into ZAP2D and would allow for the 
possibility of automatic grid manipulation (i.e., grid boundary and density) 
within ZAP2D iterations. 

A total of six grid domains have been utilized to study the blockage 
effects of the outer boundary on ARC2D and ZAP2D. A description of the 
grids and their corresponding case identifiers are listed in Table 2. In 
this table and in the remainder of this report, R is normalized by chord 
length . ° 


GRID 

R /NO. OF GRID POINTS 
o 

A 

25/185 X 65 

B 

12.5/183 X 61 

C 

5/181 X 57 

D 

1/177 X 50 

E 

0.5/177 X 45 

F 

0.25/177 X 40 


Table 2. Grid Cases. 
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A basic "C" mesh (Figure 9) around the NACA 0012 airfoil was generated by 
GRD2D which uses an algebraic method where surface normals are maintained. 
The total number of grid points range from 12,025 (185 x 65) for the largest 
domain (outer boundary characteristic distance, Ro - 25 chord lengths --Grid 
A) to 7,080 (177 x 40) for the smallest domain (outer boundary characteris- 
tic distance, Ro “ 0.25 chord lengths- -Grid F) . The upper and lower airfoil 
surface each have 80 grid points for all grids and the wake surface has ap- 
proximately 13 points. The grid 5 generated is clustered near the surface 
with a spacing of 2.5 x 10 so that at least one point can be embedded 
within a viscous sublayer. Successive grids were obtained by stripping off 
outer rings (similar to Ref. 12) to ensure that the spacial distribution 
over the airfoil remains constant for all grids. Some grid points toward 
the downstream outflow boundary have also been removed for the smaller grids 
to reduce the computational domain. Separate studies for the grids used 
here have proven that this has no effect on the calculations. Finally, the 
inner and outer boundaries for the smallest grid (F) are illustrated in 
Figure 9 . 

Numerical simulations for the NACA 0012 airfoil at two angles of attack 
are presented here: a = 4.966° and a - 10°. For both cases, the Mach num- 
ber was 0.3 and the Reynolds number was 6 million. Furthermore, the loca- 
tion of the inner boundary coupling surface, S was held fixed (at a dis- 
tance approximately 10% of the chord from the airfoil) throughout all cal- 
culations. Finally, for these calculations, the P0T2D simulation was not 
corrected for Mach number. As the local value of the Mach number at the in- 
ner and outer boundaries are quite small, little error should be introduced 
by this approximation. 


3 . 1 Validation at a «- 4.966° 

Since the inner fluid boundary surface, S. , is to approximate a poten- 
tial flow surface, it is logical to validate the coupling concept for an 
airfoil at low angle of attack and also at low Mach number. Flow separation 
and transonic effects would only force the inner boundary outward which 
would possibly diminish some of the improvement in computational effort com- 
pared to the simpler problem, but would not invalidate the coupling tech- 
nique itself. Further, the scheduling (number and frequency) of the poten- 
tial flow boundary condition updates to the viscous calculation will prob- 
ably be important as the domain size decreases. This effect should be in- 
vestigated in detail and is demonstrated in the a - 10° case, but is removed 
from the current cases illustrated in the following figures by updating only 
after ARC2D is fully converged. That is, for each domain, a total of 3,000 
iterations were utilized with ZAP2D coupling updates every 1,000 iterations. 

The effects of domain size (Grids A-F) on the computed integrated 
forces and moments are compared for ARC2D alone and ZAP2D. The ARC2D simu- 
lation is shown for both the point vortex correction (CIRCUL-TRUE) and for 
no correction (CIRCUL— FALSE) . The computational trends of ARC2D for domains 
smaller than about 4 chords (1/R^ — 0.25) are shown here for comparison only 
since these smaller domains are not recommended practice (11). Finally, 
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Figure 9. Enlarged View of NACA 0012 Airfoil Grid (185 x 65). 


ARC2D calculations are not 
code diverged for this case. 


presented for the smallest grid (F) since the 


ficiJ^ gU ? e 1 ° E . il i U ^ r f t ? S the effect of domain size on sectional lift coef- 
ncient c £ . Each division represents a change of approximately 2% in r 

compared with the experimental value as shown in th e P figure. The ZAP2D cal^ 
culatron converges to the ARC2D (CIRCUL-TRUE) solution for 1/R less than 
’ Furt J er ; 8^13 studied, C changes fro „ thls 1 ' 1?la i 

^ " 3b 1 ’ ll i UStratin 8 the of the coupled solution in 

aintainmg accuracy for grid domains down to R - 1/4 chord. The computed 
C f is also within 2% of the experimental value for all grid domains. P 

£ esult ® are illustrated in Figures 11 and 12 for sectional drag 

ter f 5wH ’r C d’ sec ^ ional Patching moment coefficient about the quar- 
ter-chord C respectively. That is, ZAP2D coupling has removed the ef- 
fects of domain size on C £ , Cj, and C m for all grids studied (A-F). 

distribution of surface pressure coefficient, C , predicted by 

FieuL^ th ® smallest 6^id (F) is compared with the experimeRtal values in 
Figure 13. The correlation is excellent, with viscous effects improving the 

and trailing edgeS ° Ver the P— tial P flow 6 case 


Contours of constant Mach number in the flow region near the airfoil 

(CIRCUL^FALSE) fo/th^l ^ ' /^ gur ? 1 ^ (a) illustrates the ARC2D calculation 
7 AP?n t V * l h largest domain (A-Grid) while 14(b) illustrates the 
ZAP2D calculation for the smallest domain (F-Grid) . The comparison between 
the two computed solutions is quite good near the airfoil with a small 
change in the ZAP2D calculation near the outer boundary of R -=0.25. 

The numerical stability and convergence characteristics are shown in 
gure where the time history of C. is compared for the first 1 000 iter- 
ations for ARC2D (large domain, A-Grd) and ZAP2D (small domain, ’ F-Grid) 

very^mall^i e" ^or^hi 0 f 0- *"**? C i, for the remaining 2,000 iterations is 
verv llttl/aFfW f ^ low . angle of attack and this grid size there is 
l °J suc ? essive Potential flow updates to the outer bound- 
ary conditions beyond setting the initial conditions. Still, the ZAP2D time 
history demonstrates that the coupling procedure is, in fact, a stabilizing 
rtHf 1 ? enCe < beCa y Se it: allows for a lfl rge reduction of grid domain. As the 
flow fielded r* information is passed more quickly throughout the 

“I - ? the oscillatory behavior evident in the ARC2D-alone simula- 

tion (Figure 15(a)) is eliminated in the ZAP2D simulation (Figure 15(b)). 

histories for the i°E the norm of the density re- 
sidual (Log(Res)) is compared for ARC2D and ZAP2D in Figure 16 The ARC2D 
large domain, A-Grid calculation is illustrated for the first MM W 
ons. urther iterations of ARC2D to a maximum of 3,000 does not improve 
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Figure 10. Computed Lift Coefficient Variation with Domain Size for a = 4.966°. 
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Figure 11. Computed Drag Coefficient Variation with Domain Size for a = 4.966°. 
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Figure 12. Computed Pitching Moment Coefficient Variation with Domain Size 
for a = 4.966°. 



Figure 13. Comparison of ZAP2D (F-Grid) and Experimental Cp Distribution 
for NACA 0012 Airfoil at Mach 0.3, a - 4.966°. 
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(a) ARC2D with A-Grid (R q = 25) 

Figure 14. Comparison of ZAP2D and ARC2D Computed Mach Contours for NACA 0012 at Mach = 0.3, a = 4.966. 
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(b) ZAP 2D with F-Grid (R q = 0.25) 
Figure 14. Concluded. 
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(b) ZAP2D with F-Grid (R q = 0.25) 


Figure 15. 
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Figure 16. Iteration History of ARC2D and ZAP2D Computed Solution Residual. 




convergence from the level indicated at 1,000. The ARC2D simulation shows a 
reduction in the residue by approximately 3 orders of magnitude. In con- 
trast, the ZAP2D small-domain, F-Grid residual is more rapidly reduced by 
almost 5 orders of magnitude due to the smaller flow domain. Abrupt changes 
in the ZAP2D calculation occur due to the outer boundary condition adjust- 
ment at the coupled iteration; however, these peaks diminish with further 
coupled iterations and the solution convergence recovers from each boundary 
update in less than 100 iterations. 

In summary, for the grids utilized in this study of the low angle-of- 
attack case has shown the following. 

1. The coupled procedure installed in ZAP2D allows for a reduction of the 
outer grid boundary from about 25 chord lengths to 0.25 chord lengths 
(R q “ 25 -» 0.25) with no loss in accuracy. 

2. The numerical convergence of the viscous ARC2D solution is greatly im- 
proved due to the reduced domain size allowable in ZAP2D. 

3. By utilizing the potential flow solution (ZAP2D) for the initial bound- 
ary conditions, a converged solution was obtained for the small domain 
in contrast to the ARC2D divergent behavior. 

4. Multiple potential/viscous boundary condition updates do not signifi- 
cantly affect the final converged solution for the smallest grid 
studied. This implies that additional grid domain size reduction 
should be possible without sacrificing accuracy. 

5. Without the further improvements in ZAP2D that may be possible such as 
accelerated convergence methods because of increased stability (Item 2 
above) or additional reduction in domain (Item 4), the ZAP2D simulation 
has still shown a reduction factor of approximately 10 in computer time 
for a given accuracy and convergence compared with ARC2D alone for this 
case. This improvement comes from the reduction in grid points and 
also a reduction of total number of iterations from 1,000 to about 200- 
300. 


3 . 2 Validation at a - 10° 

Since the a - 4.966° case did not show that multiple ZAP2D potential/ 
viscous flow updates were required to obtain the needed accuracy for the 
grids used, the a - 10° case was added to increase the extent of the viscous 
effects into the flow domain. Because of the limited time for this study, 
ARC2D was executed with CIRCUIT-FALSE only. Furthermore, only Grids D, e[ 
and F were investigated for the ZAP2D calculations since the larger grid 
cases would certainly converge to the ARC2D-alone results (Section 3.1). 
Finally, a total number of 1,000 iterations were used for all cases. 
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As early calculations indicated that intermediate updates would impact 
the convergence characteristics for this higher angle -of- attack case, a 
series of ZAP2D update schedules were investigated as outlined in Table 3. 
For example, the SI schedule includes the initial P0T2D outer boundary con- 
ditions only with no intermediate ZAP2D updates while the S6 schedule in- 
cludes nine intermediate ZAP2D outer boundary updates, 1 after each 100 
ARC2D iterations. Schedules SI through S4 were investigated with Grids D 


SCHEDULE ARC2D ITERATION NUMBER 

CORRESPONDING TO ZAP2D 
BOUNDARY UPDATE 


SI 

1 



S2 

1. 

201 


S3 

1. 

201, 

401, 601 

S4 

1. 

201, 

401, 601, 801 

S5 

1. 

101, 

201, 301 

S6 

1, 

. 101, 

, 201, 301, .. 


Table 3. Schedule of ZAP2D Outer Boundary Updates. 


and E while schedules S5 and S6 were added to the F-Grid cases when the SI 
through S4 schedules failed to converge for this small domain. This diver- 
gence of ZAP2D for the F-Grid and a - 10° case is in sharp contrast to the 
same grid and a - 4.966° where intermediate update schedules had very little 
effect on the solution. Once more, the ARC2D-alone simulation was also di- 
vergent for the smallest F-Grid, and the computational trends of the ARC2D- 
alone simulation for the smaller domains are shown here for comparison pur- 
poses only. 

The preliminary results for the effect of domain size on the sectional 
C ^ is illustrated in Figure 17. Again, each division represents a change of 
approximately 2% in C^ compared with the experimental value which is also 
shown in the figure. The ZAP2D computed C^ for Grids D and E, domain size 
reduced up to one -half chord (R =0.5), and schedule S4 deviates from the 
experimental C. and the ARc2d predicted C. for large domain (R = 25) by 
less than 1/2%. The utility of the zonal solution included in °ZAP2D in 
maintaining accuracy with reduced grid domain size is contrasted with the 
ARC2D-alone solution. On the other hand, while the ZAP2D predictions for 
Grids D and E are fully converged, the results for the smallest grid, F, are 
not completed at the present time. For this smallest grid (R = 0.25), the 
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Figure 18. Computed Drag Coefficient Variation with Domain Size for a 




outer boundary update schedule is of extreme importance at this higher angle 
of attack. As mentioned above, the ZAP2D code diverged for schedules SI 
through S4; consequently, schedules S5 and S6 were investigated and ZAP2D 
solutions for both of these later schedules are shown in Figure 17 . In this 
smaller grid domain, the effects of compressibility or the wake source term 
(see Section 2.2), which has been neglected in the potential flow updates 
thus far, may need to be included to obtain a robust ZAP2D solution. An 
early investigation of this behavior on a Cray computer also indicates that 
this problem may be related to the Unix installation of ZAP2D on the IRIS 
and TITAN workstations. Regardless, ZAP2D predicted C. for the F-Grid for 
both the S5 and S6 schedules are within approximately 5% of the experimental 
value . 

Similar results for domain dependence are illustrated in Figures 18 and 
19 for drag coefficient, C., and pitching moment coefficient about the quar- 
ter-chord, C^. The change in the ZAP2D (F-Grid) predicted C, from the ARC2D 
(A-Grid) C d is approximately 10 drag counts and the change in C compared 
with the ARC2D value is approximately 0.001 for grid domain reSuctions to 
one -half chord. On the other hand, the pitching moment predictions for both 
ARC2D and ZAP2D are approximately C - 0.01 while the experimental value is 
almost double this level. While no explanation is given for the difference 
at this time, it does illustrate that the zonal coupling in ZAP2D cannot 
improve upon a fully converged, large domain, ARC2D solution. 

The ZAP2D predicted distribution of surface pressure coefficient, C , 
for a — 10° and the E-Grid is compared with experimental measurements for ^a 
corrected c» c — 10.122“ in Figure 20. The data correlation is excellent over 
the entire airfoil with a small underprediction of the leading- edge suction 
peak. This error is most likely a result of the slightly lower angle of at- 
tack used in the numerical simulation. A similar comparison of the pressure 
distribution (not shown here) for the F-Grid indicates that the pressure 
distribution is converged except near the leading-edge suction peak where 
C (S5) - -6.25 and C (S6) - -6.38. 

"min Pm in 

A final set of three figures are included here to demonstrate the ef- 
fect of intermediate outer boundary updates (i.e., in this case, updates 
before a fully converged ARC2D solution) on the ZAP2D- E-Grid iteration his- 
tory. Figures 21 and 22 illustrate the ZAP2D solution iteration history for 
and for intermediate update schedules SI through S4. Since the SI 
schedule includes only the initial outer boundary potential flow correction, 
the impact of the various intermediate outer boundary update schedules (S2- 
S4) is evident in Figures 21 and 22. As shown, the amplitude of the C. and 
increments at the scheduled outer boundary update periods decreases ^dur- 
ing the time-step iteration history. Finally, the convergence of the inte- 
grated forces and moment for this case as a function of number of intermedi- 
ate outer boundary updates is shown in Figure 23. For all coefficients, C g , 
C^, and C^, convergence is clearly demonstrated. The percentage change In 
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Figure 19. Computed Pitching Moment Coefficient Variation with Domain Size 
for a = 10° . 


-Cp 


ZAP2D X EXPERIMENT (a = 10.12°) 


-2 


}< 

1 

\ 

y 





\ 








- 

__ 

, /v. 

"" — — 



_ - — 




x X- 

K -X-X- 

F- 

> v .--- x- 

^< xX 

--Xr — - 

x — -x — x 

X XX X 

xx : >w<Xx 


0.4 


0.2 


0.0 


- 0.2 


0.0 


0.2 


0.4 


0.6 

X 


0.8 


1.0 


- 0.4 


Figure 20. Comparison of ZAP2D (E-Grid) and Experimental Cp Distribution 
for Mach = 0.3, a = 10°. 





Figure 21. Iteration History of ZAP2D (E>*Grid) Calculated C„ (Mach = 0.3, 
a = 10°) for Various Outer Boundary Update Schedules. 



Figure 22. Iteration History of ZAP 2D (E-Grid) Calculated C (Mach = 0.3, 
a = 10°) for Various Outer Boundary Update Schedules. 
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the coefficients from their initial (SI schedule) to their final (S4 sched- 
ule) values is as follows. 

% AC * « 2 1/4% 

% AC d * 40% 

% AC - 20% 
m 

Consequently, while the percentage change in C^ is small, the change in C, 
and C^ are quite significant. 

In summary, this abbreviated study of the high angle -of- attack case has 
shown the following additional features. 

1. The coupled procedure installed in ZAP2D allows for a reduction of the 
outer grid boundary from 25 chord lengths to 0.5 chord lengths with no 
loss in accuracy. 

2. While the smallest grid, F, has not been resolved to the same accuracy 
level at this time, it has been demonstrated that additional outer 
boundary schedule updates in ZAP2D will most likely yield converged 
results within the required accuracy. Undoubtedly, the compressibility 
correction will also need to be added to ZAP2D for this case. 

3. Multiple ZAP2D potential/viscous outer boundary condition updates do 
affect the final converged solution and are required to obtain an ac- 
curate solution for all forces and moments. 

4. As the grid domain is decreased, the schedule of ZAP2D intermediate 
outer boundary updates directly affects the robustness of the ARC2D 
solution scheme (all other parameters fixed) in that a solution can be 
obtained where it was previously divergent. 
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4.0 CONCLUSIONS 


The primary objective of the work described in this report is to dem- 
onstrate the feasibility of a new potential/viscous coupling procedure for 
reducing computational effort while maintaining solution accuracy. The 
overlapping velocity coupling procedure has been used to combine a two- 
dimensional potential flow code, P0T2D, and a two-dimensional Navier-Stokes 
code, ARC2D. 

The resultant fully coupled code, ZAP2D, was utilized to compute the 
flow past a NACA 0012 airfoil at 4.966° and 10° angle of attack. The Mach 
number was 0.3 and the Reynolds number was 6 million for both cases. A 
thorough study of the numerical behavior of the computed solutions for a 
range of grid domain sizes from 25 chord lengths to 0.25 chord lengths and 
comparison with experimental data has demonstrated the following. 

1. The coupled procedure installed in ZAP2D allows for a reduction of the 
outer grid boundary from about 25 chord lengths to 0.25 chord lengths 
with no loss in accuracy. 

2. The numerical convergence of the viscous ARC2D solution is greatly im- 
proved due to the reduced domain size allowable in ZAP2D. 

3. The zonal modeling in ZAP2D increases the robustness of the ARC2D solu- 
tion scheme in that numerical divergence of ARC2D- alone is avoided for 
the smallest domain. 

4. As the Navier-Stokes computational domain is decreased and viscous (or 
Mach number) effects become more important, multiple ZAP2D potential/ 
viscous outer boundary condition updates are required to obtain ac- 
curate and converged solutions . 

5. For the a *= 4.96° case, the ZAP2D simulation has demonstrated a reduc- 
tion factor of approximately 10 in computer time for given accuracy and 
convergence compared with the ARC2D- alone simulation. With some addi- 
tional effort, the a *= 10° case should also show similar improvements. 

It is expected that by investigating additional reductions in domain 
size, the use of accelerated convergence methods because of the increased 
computational stability (Item 2 above), and optimization of AZP2D update 
schedules (Item 4 above), the computational time for the two-dimensional 
ZAP2D simulation might be reduced by a factor of approximately 20 when 
compared with the ARC2D simulation. If similar grid density reductions are 
possible in three dimensions, then the required grid for a given accuracy 
might be reduced by about 70%; i.e., only 25-30% of the large domain grid 
would be required for a fully coupled, three-dimensional solution. This 
would translate into a significant reduction in memory requirements corre- 
sponding to the reduction in number of grid points. Furthermore, additional 
improvements in the convergence behavior due to the smaller domain would 
also be possible so that perhaps one-half the number of iteration steps 
would be required with such a zonal method. Of course, additional overhead 
in computing the potential flow updates would be incurred in the three- 
dimensional simulation; however, a conservative estimate of at least one 
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order of magnitude reduction in computational time compared with the Navier- 
Stokes-only method should be achievable for a given accuracy in such a zonal 
simulation. This improvement would allow for a much more accurate three- 
dimensional simulation of the complex flow fields associated with rotary 
wing aircraft. 
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